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ABSTRACT 

According to unified schemes of Active Galactic Nuclei (AGN), the central engine is surrounded by dusty, 
optically thick clouds in a toroidal structure. We have recently developed a formalism that for the first time 
takes proper account of the clumpy nature of the AGN torus. We now provide a detailed report of our findings 
in a two-paper series. Here we present our general formalism for radiative transfer in clumpy media and 
construct its building blocks for the AGN problem — the source functions of individual dusty clouds heated 
by the AGN radiation field. We show that a fundamental difference from smooth density distributions is that 
in a clumpy medium, a large range of dust temperatures coexist at the same distance from the radiation central 
source. This distinct property explains the low dust temperatures found close to the nucleus of NGC1068 in 10 
fim interferometric observations. We find that irrespective of the overall geometry, a clumpy dust distribution 
shows only moderate variation in its spectral energy distribution, and the lOfim absorption feature is never 
deep. Furthermore, the X-ray attenuating column density is widely scattered around the column density that 
characterizes the IR emission. All of these properties are characteristic of AGN observations. The assembly of 
clouds into AGN tori and comparison with observations is presented in the companion paper. 
Subject headings: dust, extinction — galaxies: active — galaxies: Seyfert — infrared: general — quasars: 
general — radiative transfer 



1. INTRODUCTION 

Although there are numerous AGN c lasses, a unified 
scheme has been emerging s teadily (e.g., lAntonuccll 1 19931 
l2002t lUrry & Padovanilll995h . The nuclear activity is pow- 
ered by a supermassive black hole and its accretion disk, and 
this central engine is surrounded by a dusty toroidal struc- 
ture. Much of the observed diversity is simply explained as 
the result of viewing this axisymmetric geometry from dif- 
ferent angles. The torus provides anisotropic obscuration of 
the central region so that sources viewed face-on are rec- 
ognized as "type 1", those observed edge-on are "typ e 2". 
From basic considerations, iRrolik & Fje gelman (1988J) con- 
clude that the torus is likely to consist of a large number 
of individually very optically thick dusty clouds. Indeed, 
iTristram et al.l (120071) find that their VLTI interferometic ob- 
servations of the Circinus AGN "provide strong evidence for 
a clumpy or filamentary dust structure". A fundamental dif- 
ference between clumpy and smooth density distributions is 
that radiation can propagate freely between different regions 
of an optically thick medium when it is clumpy, but not 
otherwise. However, because of the difficulties in handling 
clump y media, models of the torus IR emission, beginning 
with iPier & Krolikl d 19921). utilized sm ooth density distribu- 
tions instead. |Rowan-Robinsoji| d 19951) noted the importance 
of incorporating dumpiness for realistic modeling, but did not 
offer a formalism for handling clumpy media. 

We have recently developed such a formalism and 
presented initial re ports of our modeling of AGN 
clumpy tori in INenkova. Ivezic. & E itzurl (l2002l). 
Elitzu r Nenkova. & Ivezid (l2004lh and lElitzuil d2006l 
20071) . We now offer a detailed exposition of our formalism 

1 Seneca College, Toronto, ON M2J 2X5, Canada; 
maia.nenkova@senecac.on.ca 

2 Department of Physics and Astronomy, University of Kentucky, Lexing- 
ton, KY 40506-0055; sirockmm@pa.uky.edu, moshe@pa.uky.edu 

3 Department of Astronomy, University of Washington, Seattle, WA 
98 1 05 ; ivezic @ astro.washington.edu 



and its application to AGN observations. Because of the 
wealth of details, the presentation is broken into two parts. 
In this paper we describe the dumpiness formalism and con- 
struct the source functions for the emission from individual 
optically thick dusty clouds, th e building blocks of th e AGN 
torus. In a companion paper ( INenkova et al.l 120081 part II 
hereafter) we present the assembly of these building blocks 
into a torus and the application to AGN observations. 

2. CLUMPY MEDIA 

We start by presenting a general formalism for handling 
clumpy media. For completeness, the full formalism is de- 
scribed here, including results that are utilized only in part II. 

Consider a region where the matter is concentrated in 
clouds (see Fig. [TJ. For simplicity we take all clouds to be 
identical; generalizing ou r results to a mixture of cloud p rop- 
erties is straightforward (IConwav. Elitzur. &^ Parra 2005I)- In- 
dividual clouds are characterized by their size R c , the cloud 
distribution is specified by the number of clouds per unit vol- 
ume ;ic- Denote by V c the volume of a single cloud and by <fr 
the volume filling factor of all clouds, i.e., the fraction of the 
overall volume that they occupy. We consider the medium to 
be clumpy whenever 

= « c Vc<l. (1) 

In contrast, the matter distribution is smooth, or continuous, 
when <f> ~ 1 . It is useful to introduce the number of clouds per 
unit length 

N c =n c A c = r l (2) 

where A c is the cloud cross-sectional area and I is the photon 
mean free path for travel between clouds (Fig.Q]). Since V c ~ 
A C R C , the dumpiness condition (eq.Q} is equivalent to 

(/> = A/ C fic<l, or # c <£ (3) 

The dumpiness criterion is met when the mean free path be- 
tween clouds greatly exceeds the cloud size. 
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FIG. 1 . — A region populated by clouds of size R c . Along ray s, the photon 
mean-free-path is £ and the number of clouds per unit length is Nc = t" 1 . 

2.1. Emission from a clumpy medium 

When the dumpiness criterion is met, each cloud can be 
considered a "mega-particle", a point source of intensity S c .\ 
and optical depth t\. The intensity at an arbitrary point s 
along a given path can then be calculated by applying the for- 
mal solution of radiative transfer to the clumpy medium. The 
intensity generated in segment ds 1 around a previous point 
s' along the path is S c \(s')Nc{s')ds'. Denote by Af(s' ,s) = 
f, Neds the mean number of clouds between s' and s and by 
Pesc(s',s) the probability that the radiation fr om s' will reach 
s with out absorption by intervening clouds. Natta & Panagia 
(1984) show that when the number of clouds is distributed 
around the mean J\f(s',s) according to Poisson statistics, 



P esc (s\s) = e 



where t x (s\s) = Af(s' ,s)(l -e~ Tx ); 

(4) 

the Appendix discusses the validity of the Poisson distribution 
and presents the derivation of this result. Its intuitive mean- 
ing is straightforward in two limiting cases. When t\ < 1 we 
have t\(s',s) ~ J\f(s',s)r\, the overall optical depth between 
points s' and s; that is, dumpiness is irrelevant when indi- 
vidual clouds are optically thin, and the region can be han- 
dled with the smooth-density approach. It is important to note 
that Nt\ can be large — the only requirement for this limit is 
that each cloud be optically thin. The opposite limit t\ > 1 
gives P esc (s',s) ~ e~^ (s ' s \ Even though each cloud is opti- 
cally thick, a photon can still travel between two points along 
the path if it avoids all the clouds in between. 

With this result, the intensity at s generated by clouds along 
the given ray is 



(5) 



This relation is the exact analog of the formal solution of stan- 
dard radiative transfer in continuous media, to which it reverts 
when the cloud sizes are decreased to the point that they be- 
come microscopic particles; in that case t«1 for each par- 
ticle and £~ l is the standard absorption coefficient. The only 
difference between the clumpy and continuous cases is that 
optical depth is replaced by its effective equivalent t\ and ab- 
sorption coefficient is replaced by the number of clouds per 
unit length Nc- It should be noted that, because of the partic- 
ulate composition of matter, the solution of radiative transfer 
is always valid only in a statistical sense, corresponding in 
principle to the intensity averaged along the same path over 
an ensemble of many sources with identical average proper- 
ties. However, in the smooth-density case the large number 
of dust particles along any path (a column with unity opti- 
cal depth contains ~ 10 10 particles when made up of 0.1 /im 
dust grains) implies such small relative fluctuations around 
the mean that the statistical nature can be ignored and the 
results considered deterministic. For example, the intensity 
contours produced by a smooth-density spherical distribution 



are perfectly circular for all practical purposes. In the clumpy 
case, on the other hand, the small number of particles along 
the path can lead to substantial deviations from the mean in- 
tensity /J. In contrast with the smooth density case, the in- 
tensity of a spherical cloud distribution can fluctuate signifi- 
cantly around concentric circles even though the distribution 
average properties are constant on such circles. The flux emit- 
ted from the enclosed circular area involves areal integration 
which smoothes out and reduces these fluctuations, so devi- 
ations of individual spectral energy distributions (SED) from 
their average are expected to be smaller than the brightness 
fluctuations. In general, the SED can be expected to have 
lower fluctuations than the brightness in clumpy distributions 
with a smoothly behaved function Nc- Estimating fluctua- 
tions is beyond the scope of our formalism, which is pred- 
icated on averages at the outset. An assessment of fluctua- 
tions would require an extended theoretical apparatus or anal- 
ysis of Monte-Car lo simulations, such as those described by 
iHonig et al.ld2006h . 

After characterizing the clumpy medium by the cloud size 
R c and the volume density of clouds «c we introduced an 
equivalent set of two other independent variables — the vol- 
ume filling factor <fr and the number of clouds per unit length 
Nc- Our final expression for the intensity does not involve 
<fi; only Nc enters. A complete formalism that would not in- 
voke the assumption <\> <C 1 would lead to a series expansion 
in powers of </>, and the expression we derived in eq. [5] is the 
zeroth order term in that expansion. In fact, detailed Monte 
Carlo simulations show that, to within a few percents, this ex- 
pression describes adequately clumpy media with <f> as large 
as 0.1 (Conway, Elitzur & Parra, in preparation). Since the 
intensity calculations are independent of the volume filling 
factor, their results do not provide any information on this 
quantity, nor do they provide separate information on either 
R c or nc, only on Nc- In complete analogy, the radiative trans- 
fer problem for smooth density distributions does not involve 
separately the size of the dust grains or their volume density, 
only the combination «d<Td, which determines the absorption 
coefficient and which is equivalent to Nc- 

Equation [5] accounts only for the radiation generated along 
the path by the clouds themselves. A path containing a back- 
ground source, such as the line of sight to the AGN, requires 
different handling since averaging is then meaningless. For 
such a line of sight, the intensity 4 generated when there 
are exactly k intervening clouds (k = 0, 1,2, . . .) has a Pois- 
son probability P^. The only meaningful quantities that can 
be deduced from modeling in this case are tabulations of the 
intensities and their associated probabilities P^; that is, the 
probability is Pq that the intensity is Zrj (the unobscured AGN), 
Pi that it is I\, etc. An actual source will correspond to one 
particular member of this probability distribution. 

2.2. The cloud distribution 

The only distribution required for intensity calculation is 
Nc, the number of clouds per unit length. Our interest is pri- 
marily in cloud distributions with axial symmetry in which 
Nc depends only on the distance r from the distribution cen- 
ter and the angle (3 from the equatorial plane (the comple- 
mentary of the standard polar angle). The total number of 
clouds along radial rays slanted at angle (3 is, on average, 
A/t(/3) = jNc(r,(3)dr. It is convenient to introduce as a free 
parameter the mean of the total number of clouds along radial 
rays in the equatorial plane, A/o = A/t(0). The angular profile 
A/t(/3)/A/o is expected to decline as (3 increases away from 
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the equatorial plane in AGN tori. 

2.3. Total mass in clouds 

The mass of a single cloud can be written as ~ muN c .nA c , 
where N c ji denotes column density and oth is the hydrogen 
nucleus mass. With the aid of eq. [2] and introducing the dis- 
tribution profile r]c(r,f3) = (l/Afo)Nc(r, ft) 4 , the total mass in 
clouds is 

Mc=m H N c , H Afo JrjcdV (6) 

where the integration is over the entire volume populated by 
clouds. Note that N c .nNo is the mean overall column density 
in the equatorial plane. The total mass of a smooth distribu- 
tion with the same equatorial column density would be given 
by the same expression except that the profile r/c would then 
describe the normalized gas density distribution. 

The result in equation [6] for the overall mass in the clouds 
is independent of the filling factor <\>. The only property of an 
individual cloud that enters into this expression is its column 
density jV c ,h, which is directly related to its optical depth. The 
cloud size R c is irrelevant so long as it meets the dumpiness 
criterion eq. Q 

2.4. Total number of clouds 

Similar to the intensity, the overall mass involves the num- 
ber of clouds per unit length Nc, not the cloud volume density 
«c- The latter quantity enters in the calculation of the total 
number of clouds « tot = J ncdV. Expressing « tot too in terms 
of Nc brings in the cloud size because n to t = / (Nc/A c )dV. 
With A c ~ R 2 C and R c = (j)/Nc, 

th* = J^lfidV. (7) 

Among the quantities of interest, n tot is the only one depen- 
dent upon the volume filling factor. Self-consistency of our 
formalism requires <fi <C 1 to ensure the dumpiness criterion 
and n tot 3> A/o to validate the Poisson distribution along paths 
through the source (see Appendix). Since r]c obeys the nor- 
malization jric(r,0)dr = 1, w tot is of order ~ N\/cf> 2 and the 
two requirements are mutually consistent. 

2.5. Covering Factors 

Various AGN studies have invoked the concept of a "cover- 
ing factor" in different contexts. The covering factor is gen- 
erally understood as the fraction of the sky at the AGN center 
covered by obscuring material. This is the same as the fraction 
of randomly distributed observers whose view to the center 
is blocked. Therefore the average covering factor of a ran- 
dom sample of AGN is the same as the fraction, fi, of type 2 
sources in this sample; that is, fi is the average of the fraction 
of the AGN radiation absorbed by obscuring clouds. Denote 
the spectral shape of the AGN radiation by f e \, normalized 
according to J f e \d\ = 1 . The fraction of the AGN luminos- 
ity escaping through a spherical shell of radius r centered on 
the nucleus is, on average, 

PagnO) = J d sin/3 J d\f eX P\. esc (r,P) , (8) 

where P\, esc (r, (J) is the probability for a photon of wavelength 
A emitted by the AGN in direction (3 to reach radius r (eq. |4j. 

4 The profile ?7c obeys the normalization fr/c (r, 0) dr = 1 ; note that r]c and 
Nq have dimensions of inverse length while A/o is dimensionless. 



Therefore fi = 1 - Pagn(Rom), where R oat is the torus outer 
radius, a relation that holds independent of the magnitude of 
the optical depth of individual clouds. The spectral integra- 
tion involves all wavelengths in general, and is limited to the 
relevant spectral range in cases of specific obscuring bands. 
When the clouds are optically thick to the bulk of the AGN 
radiation, P\ fisc (R ut, P) — e~^ il3 \ independent of A (eq. @). 
Utilizing the normalization of f e \ then yields 

,7r/2 

f 2 =l- e _A/T(/3) cos (3d (3. (9) 
Jo 

In a recent study iMaiolino et ai1(l2007l) defined the "dust cov- 
ering factor" as the ratio of thermal infrared emission to the 
primary AGN radiation. Since the bulk of the AGN output is 
in the optical/UV, the radiation absorbed by the torus clouds 
must be re-radiated in infrared, therefore this dust covering 
factor is the fraction fi for optical/UV. In X-rays, the frac- 
tion fi is usually derived from the statistics of sources that 
have at least one obscuring cloud, whether Compton-thick or 
not, along the line of sight to the AGN. Since the probabil- 
ity to miss all clouds is e - ^ 1 ', equation [9] holds in this case, 
too, with A/t the overall number of X-ra y obscuring clo uds, 
whether dusty or dust- free. Therefore, the Maiolinoetal] cov- 
ering factors for dust and X-rays properly correspond to the fi 
fractions for optical/UV and X-rays, respectively, and can be 
used for a meaningful comparison between similar quantities. 
Although the actual comparison is subject to many observa- 
tional uncertainties that affect differently the two wavelength 
regimes, its formal, fundamental premise is valid. 

Other definitions of a "covering factor" do not always cor- 
respond to the fraction of obscured sou rces, and sometimes 
even l ead to values that exceed unity. iKro lik & Beg elmanl 
(1988) introduce a "covering fraction" C such that the aver- 
age column density near the equatorial plane is CN c ,u- There- 
fore, their covering factor is A/o, the average number of clouds 
along radial equatorial rays. Broad and narrow line emis- 
sion from AGN is frequently calculated from an expression 
similar to eq. without considering c loud obscuration, i.e., 
in the limit t A = (e.g., lNetzer|[l990l) . The cloud distribu- 
tion is taken as spherical and the normalization of the total 
line flux is obtaine d from a "covering factor". Similar to 
IKrolik & Begelmanl this covering factor is A/ t, the total num- 
ber o f clouds along radial rays (see eq. 8 in Kaspi & Netzer 
119991) . Since a spherical distribution has fi = l -e~-^ T (eq. 
[9), this covering factor coincides with fi only when A/t -C 1 ■ 
Note that a spherical distribution with Nj = 1 would have a 
unity covering factor according to this definition, because a 
cloud is encountered in every direction on average, when in 
fact fi is only 0.63 in this case. 

The emission line covering factor is the cloud number A/t. 
not the fraction fi. This covering factor is obtained from the 
analysis of type 1 spectra under very specific approximations. 
The population of clouds dominating the calculation under 
these assumptions can be quite different from that control- 
ling the AGN obscuration, and this covering factor can differ 
a lot from the fraction fi obtained from source statistics. As 
this discussion shows, the concept of a covering factor has not 
been well defined in the literature, referring to different quan- 
tities in different contexts. When this concept is invoked, a 
proper definition requires calculation of probability from the 
detailed cloud distribution function, similar to the one for fi 
above (eq. |9}. 

3. SOURCE FUNCTION 
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FIG. 2. — Clouds at a fixed distance from the AGN at various position 
angles a with respect to the observer direction. The dust temperature varies 
on the surface of an optically thick cloud from r max on the illuminated face to 
T mm on the dark side. Therefore the emission is anisotropic and a determines 
the visible fraction of the illuminated surface area. 

The formalism just described is general. 5 Its application to 
a specific radiative process requires modeling of the source 
function for a single cloud. Here we apply this formalism to 
dusty torus clouds heated by the central AGN radiation. We 
divide the clouds into two classes. For large optical depths, 
clouds directly facing the AGN will have a higher tempera- 
ture on the illuminated face than on other surface areas (Fig. 
12). Their emission is therefore strongly anisotropic, and the 
corresponding source function S d x depends on both distance 
r from the AGN and the angle a between the directions of 
the cloud and observer 6 . Clouds whose view of the AGN is 
blocked by another cloud will be heated only indirectly by 
all other clouds. Heating of these clouds by diffuse radiation 
is, in general, much more even, and the a-dependence is ex- 
pected to be much weaker for their source function S l c x than 

for S d x . At location (r,/3), the mean number of clouds to the 
AGN is M(r, fi) = J^Ncir, fi)dr and the probability for unhin- 
dered view of the AGN is p(r,/3) = e~^ n/3 \ The general ex- 
pression for the cloud overall source function at that location 
is thus 

S ,A(r, a, 13) = p(r, p)S d cX (r, a) + [1 -p(r, /?)] S l c x (r, a, fi) (10) 

We now describe our detailed calculations of the source func- 
tions of the two classes of clouds. 

3.1. Directly illuminated clouds 

Clouds come in all shapes and forms. But whatever its 
shape, when the size of a cloud directly illuminated by the 
AGN is much smaller than r, it is indistinguishable from a flat 
patch with the same overall optical depth. Indeed, all calcula- 
tions of line emission from AGN always model the line emit- 
ting clouds as slabs illuminated by the central engine along 
the normal to their surface. However, a fundamental differ- 
ence between a flat slab and an actual cloud of the same opti- 
cal thickness is that the former presents to an observer either 
its bright or dark faces while the latter generally presents a 
combination of both. To account for this effect we construct 
"synthetic clouds" by averaging the emission from an illumi- 
nated slab over all possible slab orientations (figure [3}. This 
procedure utilizes an exact solution of radiative transfer for 
externally illuminated slabs, and we proceed now to discuss 

5 For its application to line emission from a clumpy medium see Conway 
etal2005. 

6 Denote by i the angle to the observer direction from the axis of symmetry. 
A cloud at angle fi from the equatorial plane and azimuthal angle tp from the 
axis-observer plane has 

cos a = cos /3 cos tp sin i + sin fi cos 



the properties of an irradiated dusty slab; the synthetic clouds 
are constructed from the slab solutions in 33.1.41 



O F e = L/47tr 




FIG. 3. — A synthetic cloud at distance r from the AGN and position angle 
a from the observer direction is constructed by averaging the emis sion from 
a slab with the same optical depth over all slab orientations (see eg. 1161 . The 
local illuminating flux, F c , and maximal te mperatu re on the cloud surface, 
Tmax, are uniquely related to each other (see j|3.1.2| . 



3.1.1. Slab radiative transfer 

Thanks to the planar symmetry, the slab radiative transfer 
equation can be formulated completely in optical depth space; 
neither the density profile nor the slab geometrical thickness 
are relevant. In the case of dusty slabs this invariance ex- 
tends to the temperature equation because the only attenua- 
tion of the heating radiation comes from radiative interactions. 
Therefore the slab radiative transfer problem is fully specified 
by the optical depth t\ = q\Ty, where Ty is the dust optical 
depth across the slab at visual wavelengths and q x is the rel- 
ative efficiency factor at wavelength A. The radiative transfer 
equation along a ray at angle arccos /i to the slab normal is 



A*T^ = S\(t\)-I\(t\) 
dr x 



(11) 



where t\ is measured in the normal direction from the illu- 
minated face and S\ is the source function. For dust albedo 
GJ\ and isotropic scattering, S\ = (\-TB\)B\(T) + TB\J\ where 
J\ is the angle averaged intensity and T the dust temperature, 
obtained at each point in the slab from the coupled equation 
of radiative equilibrium. 

We solve th e slab radiative transfer problem with the ID 
code DUSTY (llvezic. Nenkova. & Elitzurlll999t) . The code 
takes advantage of the scaling properties of the radiative 
transfer problem for dust absorption, emission and scattering 
(llvezic &Elitzurl 1 19971 IE97 hereafter). The solution deter- 
mines both the radiation field and the dust temperature profile 
in the slab. The dust grains are spherical wit h size distribu- 
tion from Math is. Rumpl. & Nordsieckl d!977l) . The compo- 
sition has a standard Galactic mix of 53% silicates an d 47% 
graphi te. The optical constan ts for graphite are fromlDraine 
(2003), for silicate from the Ossenkopf, Henning, & Mathis 
( 1992) "cold dust", which produces better agree ment with ob- 
serva tions of the 10 and 18 /im silicate features (ISirocky et alJ 
2008). From the optical constants DUSTY calculates the ab- 
sorption and scattering cross-sections using the Mie theory 
and replaces the grain mixture with a single-type compos- 
ite grain whose radiative constants are constructed from the 
mix average. This method reproduces adequately full calcu- 
lations of grain mixtures, especially when optical depths ar e 
large (e.g.. lEfstathiou & Rowan-Robinsodll994IWorfll2003l) . 
The handling of the dust optical properties is exact in this ap- 
proach, the only approximation is in replacing the tempera- 
tures of the different grai n component s at each point in the 
slab with a single average. Wolf (2003) finds that the temper- 
atures of different grains are within ~ ± 10% of that obtained 
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in the mean grain approximation, and that these deviations 
disappear altogether when t > 10. Henceforth, the term "dust 
temperature" implies this mean temperature of the mixture. 
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FIG. 4. — The spectral shape of th e AGN illuminating radiation for the 
standard set of parameters (see eq. |13t . 

The slab is illuminated by the AGN flux F e x at an angle f?; 
from the normal. For isotropic AGN emission with luminosity 
L, the bolometric flux at the slab position is 



(12) 
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The illuminating flux is characterized by F. and by the AGN 
normalized SED f e \ = F e \/F e . Following Rowan-Robinson 
( 1995), we employ the piecewise power-law distribution 7 



A/eA OC 




A< A h 
A h < A < A u 
A u < A < Arj 
Arj<A 



(13) 



with the following set of parameters: Ah = 0.01/im, A u = 
0.1 /im, Arj = 1/im and p = 0.5 (see fig. HJ. The effects of 
varying t he para meters from this standard set are discussed 
below in 33T31 

3.1.2. Dust temperature 

When Ty <C 1 the diffuse radiation inside the slab is neg- 
ligible. The dust temperature is the same everywhere in the 
slab and depends only on F e ; it does not depend on either Ty 
or the illumination angle 9\ (so long as the slab remains op- 
tically thin along the slanted direction). The same holds for 
any externally heated cloud as long as it is optically thin and 
its size is much smaller that the distance to the source, and 
the slab geometry faithfully mimics this aspect. By contrast, 
when the radiation source is embedded inside a smooth den- 
sity distribution (e.g., a spherical shell), the dust temperature 
does vary even in the optically thin case because of the spatial 
dilution of radiation with distance from the source. 

When ry > 1 the slab diffuse radiation contributes signif- 
icantly to dust heating. If q a \ is the efficiency factor for ab- 
sorption, let q ae = Jq a \f e \dX and q a p the corresponding av- 
erage with the Planck spectral shape. From radiative equilib- 
rium, the dust temperature at the slab illuminated boundary 
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FIG. 5. — The effects of optical thickness ry (normal to the slab sides) 
and illumination angle 9\ on the temperature of th e sl ab illuminated surface. 
Shown is the variation of the function ip (see eq. 1 14t with Ty at various 8\, 
as marked. The dust in all figures has standard Galactic composition, the 
temperature is the mixture average. The flatness of tp for 7y > 1 implies that 
the dust temperature is the same on the illuminated face of all optically thick 
slabs. 



obeys (IE97) 



aT 4 q dP (T) = -F e q ae ip(Tv,9i) 



(14) 



The function ip introduced here accounts for the contribution 
of the diffuse radiation to the energy density, and can be ob- 
tained only from a detailed solution of the radiative transfer 
problem; it is normalized according to ip(0, 6>;) = 1 to recover 
the temperature of optically thin dust (IE97 8 ). Irrespective 
of the actual value of the dust temperature on the illuminated 
face, the function ip fully describes the effects of illumina- 
tion angle and slab overall optical depth on that temperature. 
Figure [5] shows the variation of tp with both Ty and #j. The 
dependence on either variable is weak. Even for normal il- 
lumination, where it is the largest, tp does not exceed 1.35, 
reaching its maximum at ry > 1; i.e., the contribution of the 
diffuse component to the surface heating is no more than 35% 
of the direct heating by the external source. This behavior is 
markedly different from the spherical case where ip increases 
monotonically with r v without bound (IE97). The reason is 
the fundamental difference in the relation between radiative 
flux and energy density in the two geometries. In a spherical 
shell, the energy density can increase indefinitely in the inner 
cavity because there is no energy transport there (from sym- 
metry, the diffuse flux vanishes in the cavity). As the optical 
depth increases, the radiation trapped in the cavity dominates 
the dust heating on the shell inner boundary, leading to an 
unbounded increase of tp with ry (the same surface temper- 
ature requires a smaller F e ). In contrast, a slab cannot store 
energy anywhere. When a slab optical thickness increases, 
any increase in diffuse energy density is accompanied by a 
commensurate flux increase, producing the saturation effect 
evident in the figure. As a result, the temperature on the illu- 
minated face is the same, independent of optical depth when 

Figure[6]displays the temperature variation inside represen- 
tative slabs of varying optical thickness heated to a surface 
temperature of 850 K. After a sharp drop near the illuminated 
face, the temperature is nearly constant throughout the rest 



Note that F v oc v a implies XF X oc A~ (a+1) 



Note that IE97 utilized the function <J/ = (q^/q^ip 



6 



Nenkova et al 



1000 

800 

600 
I 

400 A 

Vs. 

200 
1000 



0.1 



a = o° 



io__ 

ST— 4.fi~ 

500365= 




9 : =80 L 



10__ 

lOfC^OfTTdocf 



0.0 



0.2 



0.4 



T / T 



0.6 
T 



O.i 



1.0 



FIG. 6. — Temperature variation inside slabs of various ry, as marked, for 
two illumination angles. The variable t/t t is thickness into the slab from the 
illuminated face, where the temperature is 850 K for all the displayed slabs. 
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FIG. 7. — An optically thick dusty slab at distance r from an AGN with 
luminosity L45 = L/(10 45 erg s -1 ) is irradiated normal to the surface. The dust 
temperature is T max on the slab illuminated face, T m \ n on the dark side. The 
curve for r max is applicable whenever ry > 1, the one for 7" ra j n when ry > 10. 

of the slab for optical depths exceeding tv ~ 10. There is 
little difference between the profiles at large overall optical 
depths — all slabs with tv > 100 produce nearly indistinguish- 
able profiles. Increasing 6{ has a similar effect to increasing tv 
because of the slanted illumination. The temperature decline 
near the illuminated surface reflects the exponential attenua- 
tion of the external heating radiation. Once this radiation is 
extinguished, the dust temperature in the rest of the slab is 
maintained by the diffuse radiation. And since the bolometric 
flux is constant and must be radiated from the slab dark side, 



this temperature is uniquely determined by the external flux 
and the grain optical properties. It is independent of the over- 
all optical depth once that optical depth exceeds unity for all 
wavelengths around the Planck peak. This behavior is in stark 
contrast with all centrally heated dust distributions where the 
flux radiated from the outer boundary decreases as the surface 
area increases with overall size, leading to ever declining dust 
temperatures. 

As is evident from figure [5] ip(Ty,0i) is nearly constant as 
a function of tv for Ty>l; deviations are within the 5% ac- 
curacy of our numerical results. Therefore, at a given dis- 
tance r from the AGN and a given inclination 6 U all optically 
thick slabs are heated to the same surface temperature. And 
as is evident from figure[6] they also have essentially the same 
dust temperature on their dark face when ry > 10. Figure [7] 
displays the variation of these two surface temperatures with 
distance from the AGN for optically thick slabs with normal 
illumination. The displayed relations can be described with 
simple analytic approximations. The temperature r max on the 

1 /2 

illuminated face reaches 1500 K at distance R& ~ 0.4L 4 5 pc, 
where L45 = L/(10 45 erg s" 1 ). Introduce p = r/R^, then 

f 1500-(l/p) 039 K 



■ < 1 



640-(9/p) u - 45 K p>9 



400-(l/p) 042 K 



(15) 



These temperatures bracket the range of dust temperatures on 
the surface of a cloud at distance p from the AGN (see figure 
|2]). Because the dust temperature is determined by radiative 
equilibrium, the large disparity between different surface ar- 
eas cannot be equalized by cloud rotation. The dust is pri- 
marily heated by absorption of short wavelength radiation in 
electronic transitions, with typical time scales of 10~ 6 sec, and 
cools via vibrational transitions within seconds. Therefore 
dust radiative equilibrium is established instantaneously when 
compared with any dynamical time scale. It is also important 
to note that the gas has no effect on the dust temperature as 
long as the density is < 10 12 cm" 3 . 

All the clouds illustrated in figure [2] are at the same dis- 
tance from the AGN, yet their radiation toward the observer 
involves a mix of a wide range of surface temperatures. This 
result has profound implications for the torus emission. The 
10 pm dust emission in N GC1068 has been r ecently resolved 
in VLTI interferometry by iJaffe et all (12004), who analyzed 
their observations with a model containing a compact (r < 
0.5 pc) hot (> 800 K) core and coo ler (320 K) dust extending 
to r ~ 1 .7 pc. iPoncelet et al.l ((2006) reanalyzed the same data 
with slightly different assumptions and reached similar con- 
clusions — the coolest component in their model has an aver- 
age temperature of 226 K and extends to r = 2.7 pc. As noted 
by the latter authors, the close proximity of regions with dust 
temperature of > 800 K on one hand and ^200-300 K on the 
other is a most puzzling, fundamental problem. Clumpiness 
provides a natural solution, thanks to the spatial collocation of 
widely different dust temperatures. With a bolometric lumi - 
nosity of - 2 x 10 45 erg s" 1 for NGC1068 dMason et alj2006l) . 
the dust temperature in an optically thick cloud at r = 2 pc is 
960 K on the bright side but only 247 K on the dark side, 
declining further to 209 K at 3 pc. Indeed, the temperatures 
deduced in the model synthesis of the VLTI data fall inside 
the range cover ed by the cloud surface te mperatures at the de- 
rived distances. Schartm ann et alj d2005l) have recently mod- 
eled the NGC1068 torus with multi-grain dust in a smooth 
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density distribution and found that, although the different dust 
components span at the same location a range of temperatures, 
this range is smaller than required. They conclude that even 
with this refinement their model cannot explain the VLTI in- 
terferometric observations and that the clumpy structure of 
the dust configuration must be included in realistic model- 
ing. The same effec t is found in the Circinus AGN, where 
Tristra rrTet al.l ([2007) conclude similarly that the temperature 
distribution inferred from their VLTI observations provides 
strong evidence for a clumpy or filamentary dust structure. 

At a given distance r, the highest surface temperature is 
obtained for slab orientation normal to the radius vector. This 
maximal temperature has a one-to-one correspondence with 
F e . Thanks to this correspondence, T max can be used instead 
of F e to characterize the boundary condition of the solution. 

3.1.3. Emerging intensity from a slab 

The spectral shape of the radiation emerging from the slab 
depends on ry, r max , 0{ and the angle 9 between the observer 
direction and the slab normal (see figure[3]l. Figures l8lfToldis- 
play the dependence of the slab SED on these four parameters. 
Each figure displays separately the SEDs for the illuminated 
and dark sides of the slab, which are fundamentally different 
from each other. Every SED is scaled with F e , the AGN bolo- 
metric flux at the slab location (see eq.[T2l. 
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FIG. 8. — Optical depth dependence of SED of slabs illuminated by normal 
radiation to r raax = 850 K. In the top and bottom panels the observer direction 
is 60° from slab normal on the illuminated and dark sides, respectively. 

Figure [H] displays the Ty dependence of the SED for rep- 
resentative parameters. On the illuminated side the SED is 
dominated by scattering and hot emission from the surface 
layer. The silicate 10/im feature appears always in emission, 
although its contrast decreases with Ty. On the dark side, the 
feature displays the behavior familiar from spherical models 



(e.g., IE97): emission at small ry, switching to absorption at 
large ry when the hot radiation from the illuminated region 
propagates through optically thick cool layers. However, the 
absorption feature is never as deep as in the spherical case , 
reflecting the flat temperature profile (fig. [6]). 
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FIG. 9. — The SED of slabs with ry = 100 normally illuminated to surface 
temperature r max , as marked. All other parameters are the same as in Fig. [8] 

Figure [9] shows the variation of the SED with the temper- 
ature of the illuminated surface. The change in temperature 
affects the dust emission, shifting it to longer wavelengths as 
T decreases. This modifies the 10/im silicate feature. When 
the emission peak moves past 10/im, the silicate feature starts 
disappearing from the SED. Short wavelengths (A < 1 /im) are 
dominated by the scattered component. They are visible only 
on the bright side and are unaffected by this change. 

Figure[l0]shows the effects of the illumination and viewing 
angles. Moving the observer's direction away from the slab 
normal has a similar effect to increasing the slab optical depth. 
Varying the illumination angle affects the attenuation of the 
external heating radiation. 

3.1.4. Cloud source function x 

A slab-like patch observed at angle 9 Q from a large distance 
appears as a point whose intensity (flux per solid angle) is 
h(9o)cos6 . At a distance r from the AGN, corresponding 
to temperature r max , and position angle a we construct the 
source function for a "synthetic cloud" with optical depth ry 
from 

S^ x (T m!ix ,TY,a)= (I(T m!i x,Tv ,8i,9 )cos8 ) (16) 

Here I\(T msLX ,Ty,6i,9 ) is the intensity of a slab with the listed 
parameters and the brackets denote averaging over all possi- 
ble slab orientations (Fig. [3]). The fraction of slabs with ob- 
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FIG. 10. — Dependence of slab SED on illumination and observation an- 
gles 9\ and 9 a . The observer is on the illuminated side in the left panels, 
on the dark side in the right panels. All slabs have ry = 100, the external 
flux F c corresponds to r m ax = 850 K (the surface temperature that a normally 
illuminated slab would have at that location). 
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FIG. 11. — Dependence of the source function of directly illuminated 
clouds on optical depth. The clouds are located at r-equivalent temperature 
Tmax = 850 K and position angles as marked; a = 0° is directly in front of the 
AGN, a = 180° directly behind (see fig.0. 

servable bright side is the same as the visible fraction of the 
illuminated area on the surface of a spherical-like cloud. 

Figure QT| displays the dependence of the source function 
5^ A on optical depth for clouds located at the same distance 
at a number of representative position angles around the AGN. 
At a = 0°, only the dark side of the cloud is visible (see fig. 
|2]). Increasing a exposes to the observer a larger fraction of 
the illuminated area until it is fully visible at a = 180°. This 
explains the emergence of short wavelengths and the switch 
from absorption to emission of the 10 silicate feature as 
a increases. The same behavior is evident also in figure [121 
which provides a more detailed coverage of the a-dependence 
of Sjf A as well as additional temperatures. 

At all position angles, increasing the optical depth of a sin- 
gle cloud beyond tv ~ 100 has only a minimal effect on the 
spectral shape. And except for the short wavelengths at a = 0, 
the SED similarity extends down almost to clouds with tv 
= 10. This behavior reflects the saturation of the slab tem- 
perature profile (see fig. |6]l and is a realistic depiction of the 
situation inside an actual cloud heated from outside by a dis- 
tant source. An important consequence of the flatness of the 
temperature profile through most of the cloud is that the sili- 



cate absorption feature never becomes deep — a uniform tem- 
perature source cannot produce an absorption feature at all. 
This is fundamentally different from continuous dust distri- 
butions surrounding the heating source, where the tempera- 
ture keeps decreasing with radial distance and the absorption 
feature keeps getting deeper as the overall optical depth in- 
creases until finally the entire wavelength region is suppressed 
by self absorption. In an externally heated cloud, on the other 
hand, the depth of the absorption feature is set by the con- 
trast between temperatures on the two faces. When both T max 
and r m i n are higher than ~ 300 K (the Planck-equivalent of 
10/im, roughly) the entire slab contributes to 10/im emission; 
when both are lower, no region emits appreciably at this wave- 
length. The absorption feature is deepest when the hot face is 
warmer than ~ 300 K and the cool side is cooler, maximizing 
the contrast between the emitting and absorbing layers. The 
solutions displayed in figure [12] for T mdx = 400K and 200K 
show the deepest absorption features produced by directly il- 
luminated clouds for any combination of the parameters. 

3.1.5. The input spectral shape 

All previous descriptions of the spectral shape of the AGN 
input continu um are summarized by the piecewise power law 
in eq. [T31(see lRowan- Robinso nl 1995L and references therein). 
The continuum shape shortward of 0.1 /im (1000 A) is uncer- 
tain; the slope of the spectral falloff between optical and X- 
rays has been determined for many sources but the location 
Ah of the turndown from flat Xf e \ toward high frequencies re- 
mains unknown. Fortunately, this spectral region has no effect 
on the shape of the dust SED. It affects only the normalization 
of the relation between luminosity and dust temperature (eq. 
[15] ), and we find this to be only ~ 1 % effect when Ah increases 
from our standard 0.01/im all the way to 0.03/mi. 

The impact of the optical — near-IR region is much more 
significant. Since dust cannot be hotter than its sublimation 
temperature (~ 1500 K), it emits predominantly at A > 2- 
3/im. Shorter wavelengths involve the Wien tail of emission 
by the hottest dust and scattering of the AGN radiation, and 
thus reflect directly the input spectrum. We characterize the 
AGN emission in this region by the spectral index p and the 
wavelength Art, which marks the onset of the Rayleigh-Jeans 
tail / e A oc A -4 . This wavelength corresponds roughly to the 
lowest temperature on the accretion disk of the central engine. 
The spectral index p has been determine d for 6868 quasar s 
studied in the Sloan Digital Sky Survey (llvezic et al.1 12002). 
Its distribution covers the range -1 < p < 1.5 and has a flat 
peak at 0.5 < p < 0.8. Direct observational determination of 
Art requires near-IR spectral studies of type 1 AGN with an- 
gular resolution better than 0."01 for a source distance of 10 
Mpc. Such observations are yet to be performed. Promising 
indirect methods for separating the near-IR emission into its 
torus and disk co mponents include mea surements of contin- 
uum polarization (Kishi moto et alj2005[) and multip le regres- 
sion analysis of time variability dTomita et alj|2006l) . 

Figure [T2l shows SEDs for clouds illuminated by radiation 
with p = 0.5 and p = 1, two cases that bracket the peak of the 
observed p-distribution. The input spectrum makes no differ- 
ence at wavelengths dominated by dust emission, A > 2pm. 
Shorter wavelengths are dominated by the AGN scattered ra- 
diation, and the displayed results reflect the difference in input 
radiation in this spectral region. Although not large, these dif- 
ferences are important because IR spectral indices that include 
at least one wavelength sh orter than 2/xm are a comm on tool 
of spectral analysis (e.g.. lAlonso-Herrero et alj |2003). The 
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figure also shows an intermediate slope, p = 0.75, with Arj 
increased from 1/im to 4pm. As is evident from the figure, 
such an increase will also produce an observable effect on IR 
spectral indices. 
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FIG. 13. — Temperature variation inside slabs of various Ty, as marked, 
heated indirectly by diffuse radiation at distance from the AGN with equiva- 
lent temperature r max = 850 K (cf fig. [6). The variable t/t t is thickness into 
the slab from one face. 



3.2. Clouds Heated Indirectly 

Clouds whose line-of-sight to the AGN is blocked by an- 
other cloud will be heated only indirectly by the diffuse ra- 
diation from all other clouds. Just as in the standard, smooth 
density case, the self-consistent solution for the coupled prob- 
lems of diffuse radiation field and source function can be read- 
ily obtained by A-iterations. In the first step, calculate the 
source functions for all directly heated clouds, and from the 
emission of these clouds devise a first approximation for the 
diffuse radiation field. Next, place clouds in this radiation 




100 



FIG. 14. — Source functions of clouds heated indirectly by the diffuse radi- 
ation for two optical depths and a set of distances, labelled by the equivalent 
temperature r max . 

field and calculate their emission to derive a first approxima- 
tion for the source functions of indirectly illuminated clouds 
and, from eq.[10] the composite source function at every loca- 
tion. In successive iterations, add to the AGN direct field the 
cloud radiation calculated from eq. [3] and repeat the process 
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until convergence. The solution must be tested for flux con- 
servation, which takes the same form in both the smooth and 
clumpy cases: JdXj F It \(r,£l)d£l/4Tv = F e , where the angu- 
lar integration is over a spherical Gaussian surface of radius 
r and F r ,\ is the radial component of the radiative flux vector, 
comprised of the AGN transmitted radiation and the diffuse 
dust emission. In clumpy distributions, the contribution of the 
transmitted flux to this integral becomes /?AGN^e, where /?agn 
is from eq. [8] and the radial component of the diffuse flux at 
any position is F^ x = J I^cosddft, where l9 is the intensity 
from eq. [5] and 9 is angle to the radius vector. Thus the flux 
conservation relation becomes 



10' 



PAGN^e + J d sin/3 Jd\F r c x (r, f3) = F e 



(17) 



With trivial modifications, this is also the result for the 
smooth-density case: in the transmitted term (pagnF^) the 
clumpy effective optical depth t\ is replaced by actual opti- 
cal depth t\ in the escape factor Px, esc , and in the diffuse term 
7 A stands for the smooth-density emission intensity. 

Dust is heated predominantly at short wavelengths and is 
reprocessing radiation toward longer wavelengths. As a re- 
sult, at every location the clouds heated directly by the AGN 
will be the warmest and provide much stronger heating than 
the shadowed, much cooler clouds. The significance of feed- 
back from the clouds that are heated indirectly can be gauged 
from S l c A /SJ? A , the ratio of their source function to the driving 

term Sj! A ; when this ratio is small, rapid convergence can be 
expected. In our calculations we performed only the first two 
steps of the A-iteration process and employed an isotropic ra- 
diation field as an approximation for the heating of shadowed 
clouds. The approximate heating field was derived from an 
average over a of the emission of directly heated clouds at the 
given location; this is the radiation field that would exist inside 
a spherical shell of directly illuminated clouds at the given lo- 
cation. In reality, an indirectly heated cloud is probably ex- 
posed to the bright sides of somewhat fewer clouds, because 
they are on the far side from the AGN, so our approximation 
is an upper limit for the strength of the diffuse field. In this 
isotropic field we placed dusty slabs and solved the radiative 
transfer problem with DUSTY. Figure D~3lshows the tempera- 
ture profile inside such slabs at the location where T mlLX = 850 
K. As expected, slabs heated by the diffuse radiation are much 
cooler than their directly illuminated counterparts, shown in 
fig. [6] As before, from the radiative transfer solution we ob- 
tained S l c A using the averaging procedure in eq.[l6] Since now 
the slab is illuminated isotropically on both faces, the depen- 
dence on 9i disappears and S l c A is isotropic. Figure FBI shows 

results for SJ, A at two representative optical depths and a range 
of radial distances, characterized by T m . dx . Figure Q3] shows 
the corresponding ratio S l c A /5^ A for tv = 100. This ratio is 
always below ~10% around the wavelengths corresponding 
to peak emission at each distance; for example, clouds with 
?max = 200 are the main contributors to the emission around 
~15^m, and as is evident from the figure, Sj. A is less than 

10% of x for this temperature in that spectral range. Sig- 
nificantly, S l c A /Sj? A is small even though our approximation 
overestimates the diffuse radiation strength as it involves an 
isotropic distribution of directly heated clouds around every 
AGN-obscured cloud. As expected, clouds heated indirectly 
are much weaker emitters, indicating rapid convergence of the 
iteration procedure. 
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FIG. 1 5 . — Ratio of the source functions for directly- and indirectly-heated 
clouds with ry = 100 at the indicated distance-equivalent temperatures. S A is 
taken at a = 90° . 

In contrast with figs. QT| and [12] figure [14] shows only a 
trace of the 10 /im feature. This disparity reflects the fun- 
damental differences, evident in figures [6] and Q~3] between 
the internal temperature structures of directly- and indirectly- 
illuminated clouds. A cloud placed in a black-body radia- 
tion field will thermalize with its temperature T and emit ac- 
cording to I\ = B\(T)(\ -e~ Tx ). Therefore, when rio^m <S 1, 
I\ ~ B\T\ and the emergent spectrum has the same shape as 
the dust cross section, producing an emission feature. When 
Tio/xm increases and approaches unity, self-absorption sets in 
and the feature strength decreases. As Tio/im increases further 
and exceeds unity, I\ becomes equal to B\\ at constant tem- 
perature, self-absorption and emission exactly balance each 
other, producing the thermodynamic limit of the Planck func- 
tion. For all Tio^m > 1, a single-temperature cloud will never 
produce a feature, either in emission or absorption. Since our 
diffuse radiation field differs from a pure black-body, the tem- 
perature inside indirectly heated clouds is not constant, but 
its variation is still relatively modest, as is evident from fig. 
n~3l While the dust temperature in a directly illuminated slabs 
varies by more than factor 4, it varies by less than factor 2 
inside slabs heated from both sides by the isotropic diffuse ra- 
diation field. Furthermore, even this small variation is mostly 
limited to narrow regions near the heated surfaces. This ex- 
plains why fig. [T4l shows only a weak emission feature at tv 
= 10 (which corresponds to Tio^m — 0.7) and practically no 
feature in most of the tv = 100 clouds. 

4. DISCUSSION 

The formalism presented here is general and can be applied 
to any clumpy distribution. In this study we employed it to 
construct the source functions of dusty clouds. As a small ob- 
ject, a cloud is primarily characterized by its overall optical 
depth tv. Two additional properties could potentially affect 
the cloud emission. One is surface smoothness. A rough, 
fractal-like surface can be expected to reduce the efficiency 
of scattering by absorbing some scattered photons and trans- 
ferring them to the thermal pool. The other is shape. This 
factor could be studies by averaging ellipsoidal clouds over 
all orientations. The shape parameter would then correspond 
to the ellipsoid axial ratio as the cloud varied from a spherical 
shape to extreme elongation; the clouds constructed here can 
be considered representative of extremely elongated shapes. 
The isotropy of the external radiation field of indirectly heated 
clouds enabled us to explore in this case the two extremes 
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of cloud geometrical shape: since the radiative transfer prob- 
lem for a spherical cloud retains the spherical symmetry in an 
isotropic external field, it too can be solved by DUSTY. The 
source function is then found from SJ. A = F\/Q, where F\ 
and O are, respectively, the flux and solid angle of the sphere 
at a large distance. We calculated S]. A for spherical clouds 
imbedded in the same radiation field used in the calculations 
described in 33.21 The solution for a sphere with uniform 
density depends only on the overall optical depth ry, allow- 
ing direct comparison with the clouds constructed by averag- 
ing slabs of the same ry. We found no significant differences 
between the two extreme cases. We are currently studying 
also the direct illumination of spherical clouds, and will report 
our full findings in a separate publication (Kimball, Ivezic & 
Elitzur, in preparation). 

Both shape and surface properties appear to be of secondary 
importance. The modeling of the source function presented 
here in which the clouds are characterized purely by Ty does 
seem to capture the essence of emission from a single dusty 
cloud, although definitive conclusions will have to await de- 
tailed comparisons with the calculations for externally illumi- 
nated spheres. 

4.1. Dust temperature in clumpy media 

Dust temperature distributions are profoundly different in 
clumpy and smooth media. In smooth density distributions, 
dust temperature and distance from the AGN are uniquely re- 
lated to each other — given the distance, the dust temperature 
is known, and vice versa. In contrast, as shown in 33.1.21 in a 
clumpy medium 

• different dust temperatures coexist at the same distance 
from the AGN 

• the same dust temperature occurs at different distances 
— the dark side of a cloud close to the AGN can be as 
warm as the bright side of a farther cloud 

ISchartmann et al.l (|2005) find that different components of the 
multi-grain mix employed in their torus modeling can have 
different temperatures at the same location, therefore the con- 
cept of "dust temperature" is ill defined at the microscopic 
level in such regions. Still, the fundamental differences be- 
tween the temperature structures of smooth and clumpy media 
apply to the individual components of the grain mix and to the 
mixture average, as well as to the common, proper dust tem- 
perature when all the components equilibrate to the same one. 
These differences have profound implications for the torus 
emission, explaining the low dust temperatures found so close 
to the nuclei of NG C1068 (iJaffe et alJl2004t iPoncelet et all 
120061) and Circinus dTristram et al.ll2007l) . We discuss these 
implications further in part II. 

4.2. X-rays vs IR 

IR flux measurements collect the emission from the entire 
torus area on the plane of the sky. This flux is determined 
by the average number of clouds along all radial rays through 
the torus. In contrast, X-ray attenuation is controlled by the 
clouds along just one particular ray, the line of sight to the 
AGN. Since X-rays are absorbed by dust-free as well as dusty 
material, X-ray absorbing clouds will generally outnumber 
the torus clouds in any given direction. In fact, X-ray ab- 
sorption in AGN could be dominated by the dust-free clouds 
(see paper II). But even for the dusty portion of the column, 



the number of X-ray absorbing clouds can differ substantially 
from the torus average. As an example, the appendix table 
lAl presents a tabulation of the Poisson distribution for Af = 
5, which is representative of AGN tori as shown in part II. 
More than 80% of the paths will have a number of clouds dif- 
ferent from 5 in this case, and the probability to encounter 
just 1 cloud or as many as 9 is a full 20% of the probability 
to encounter the average 5. Two type 2 sources with similar 
cloud properties and the same average Af will have an iden- 
tical IR appearance, yet the X-ray absorbing columns in each 
torus could still differ by an order of magnitude. This can be 
expected to introduce a large scatter in X-ray observations. 

Each spectral regime responds to large variations in cloud 
optical depth in an entirely different way. The IR emission 
depends on Ty through the probability for photon escape and 
the cloud source function (see eq. 0. Both factors saturate 
when Ty exceeds ^100. From eq.[4] f e sc = <? _Af when t\ ^> 1. 
Since this condition is met at all relevant wavelengths when 
ry > 50, Pesc becomes independent of ry. Similarly for S C X- 
Because each cloud is heated from outside, only its surface is 
heated significantly when ry is large. Increasing ry further 
only adds cool material, thus S c .\ saturates for all relevant A 
(similar to standard black-body emission). Indeed, figure [TT] 
displays model results that cover three orders of magnitude of 
clump optical depth, yet the SEDs show only moderate varia- 
tion, saturating altogether when ry > 100. Even at ry < 100, 
significant spectral variety is mostly limited to A < 10 /im 
for clouds along the line of sight to the AGN. In contrast, 
the Ty-dependence of X-ray absorption is markedly different. 
Individual torus clouds are optically thin to X-rays, because 
the optical depth for Thomson scattering is only ~2x 10~ 3 Ty, 
therefore the overall optical depth for X-ray absorption is 
A/Yy (see 32.11 ). It increases linearly with ry, in contrast with 
the saturated response in the IR regime. 

The great disparity between the two spectral regions is ex- 
pected to further increase the scatter in torus X-ray proper- 
ties among AGN with similar IR emission. It also may help 
explain why the SEDs show only moderate variations in the 
infrared that a re not well correl ated with the X-ray absorbing 
columns (e.g.. lSilva et al.ll2004l) . 

4.3. The 10\im feature 

In contrast with ultra-luminous infrared galaxies (ULIRGs), 
AGN observations do not pr ovide any examp le of deep 10/im 
silicate absorption feature (lHao et al.l 120071) . This behav- 
ior conflicts with th e results of smooth-density models (e.g., 
iPier & Krolikl [l992) but is a natural consequence of dumpi- 
ness: as is evident from figuresQT|and[T2l single clumps never 
produce extremely deep features. This behavior reflects the 
flat slab temperature profile (see fig. |6]) and is a realistic de- 
piction of the situation inside an actua l cloud heated from out - 
side by a distant source. As noted by iLevenson et al.l (2007), 
the different behavior of the 10/im absorption in ULIRGs and 
AGN indicates that the dust distributi on is smooth in th e for- 
mer and clumpy in the latter (see also Spoon et al. 200~7|). 

Clumpiness suffices by itself to explain the modest depth of 
the 10 /im absorption feature in AGN. The complete behavior 
of the feature, including the transition from emission to ab- 
sorption, involves an intricate interplay between the relative 
contributions of clouds at different locations and their mutual 
shadowing. This behavior displays a complex pattern that de- 
pends on the actual geometry of the cloud distribution. A de- 
tailed, quantitative analysis of the 10/im feature in clumpy tori 
is performed in part II. 
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4.4. Conclusions 

Clumpy media differ from continuous ones in a number of 
fundamental ways. The low dust temperatures found close to 
the nucleus of NGC1068 contradict basic physical principles 
in smooth density distributions but arise naturally in clumpy 
ones. Two additional puzzling features of IR emission from 
AGN are simply hallmarks of clumpy dust distributions, inde- 
pendent of the distribution geometry: Even if the cloud con- 
figurations were spherical or irregular rather than toroidal, (1) 
the SED still would show only a moderate range incommensu- 
rate with the variation of X-ray attenuation, and (2) the 10/im 
absorption feature would never be deep. Understanding the 



full range of 10/im features observed in AGN spectra requires 
considerations of the actual geometry of the cloud distribu- 
tion, though. This problem is addressed in part II together 
with other implications of clumpy tori to AGN observations. 



Part of this work was performed while M.E. spent a most 
enjoyable sabbatical at LAOG, Grenoble. We thank Nancy 
Levenson, Paulina Lira and Hagai Netzer for their useful com- 
ments on the manuscript. Partial support by NSF and NASA 
is gratefully acknowledged. 



APPENDIX 
POISSON STATISTICS 



TABLE 1 



* Pk Pk/PAf T, p i 
i<k 






0.0067 


0.04 




1 


0.0337 


0.19 


0.04 


2 


0.0842 


0.48 


0.12 


3 


0.1404 


0.80 


0.27 


4 


0.1755 


1.00 


0.44 


5 


0.1755 


1.00 


0.62 


6 


0.1462 


0.83 


0.76 


7 


0.1044 


0.60 


0.87 


8 


0.0653 


0.37 


0.93 


9 


0.0363 


0.21 


0.97 


to 


0.0181 


0.10 


0.99 


11 


0.0082 


0.05 


0.99 


12 


0.0034 


0.02 


1.00 



NOTE. — Poisson probability (second column) for the number k listed in the first column when the average is Af = 5. The third column normalizes these 
probabilities to the one for hitting the average Af. The last column lists the cumulative probability. 



The elementary problem most relevant for the statistics of clouds along a ray is the distribution of points placed randomly on 
a circular board. Denote by n the total number of points and by p (<C 1) the probability that a point lands on a given radial ray. 
Then the probability for k points to land on that ray is 



k\(n-k)\' 

If Af is the average number of points landing on a ray than p = Af/n, therefore 



Af k 



n 



1-1 
n 



I- 2 - 

n 



k+l 



n 



(Al) 



(A2) 



e ^ , yielding the 



In the limit n — > oo with Af and k fixed, every term in the square brackets approaches unity while (1 -Af/n)" 
Poisson distribution 

* - kl • 

Equation 2] for the photon escape probability follows immediately. The probability to encounter k clouds along the path is If 
the optical depth of each cloud is r then the probability to escape from all of these encounterers is e~ kT . Therefore 



(A3) 



k 



(Afe~ T ) k 



(A4) 



yielding the Natta & Panagia result (eq. |4]l- 

It is important to note that the only requirement for eq. |A3| is that the total number of points obey n k,Af; the mean number 
Af of points on a ray can be small (even less than unity). While the statistics of points along a ray always follows the Poisson 
distribution, there is no similar universal limit for the statistics of the average Af, and it remains arbitrary. As an example, consider 
a large number of identical boards with the same number of points n spread on each one of them, so that the average Af is the 
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same for every board. In this case TV will have a <5-function distribution while the number of points on any given ray in each 
board will obey the Poisson distribution around that common average. 

Finally, the Poisson distribution allows significant deviations from the average. From the accompanying tabulation for J\f = 
5, the probability to hit just one cloud or as many as nine clouds is 20% of the probability to hit the average five clouds. This 
implies a substantial probability that the number of clouds along any particular line of sight, such as the one to the AGN, will 
deviate significantly from the torus average /V. 

REFERENCES 



Alonso-Herrero, A., Quillen, A. C, Rieke, G. H., Ivanov, V. D., & Efstathiou, 

A. 2003, AJ, 126, 81 
Antonucci, R. 1993, ARA&A, 31, 473 

Antonucci, R. 2002, in Astrophysical Spectropolarimetry, ed. J. Trujillo- 

Bueno, F. Moreno-Insertis, & F. Sanchez, 151-175 
Conway, J., Elitzur, M., & Parra, R. 2005, Ap&SS, 295, 319 
Draine, B. T. 2003, ApJ, 598, 1017 

Efstathiou, A., & Rowan-Robinson, M. 1994, MNRAS, 266, 212 
Elitzur, M. 2006, New Astronomy Review, 50, 728 

Elitzur, M. 2007, in ASP Conf. Ser. 373: The Central Engine of Active 

Galactic Nuclei, ed. L. C. Ho & J.-M. Wang, 415^124 
Elitzur, M., Nenkova, M., & Ivezic, Z. 2004, in ASP Conf. Ser. 320: 

The Neutral ISM in Starburst Galaxies, ed. S. Aalto, S. Huttemeister, & 

A. Pedlar, 242-252 
Hao, L., Weedman, D. W., Spoon, H. W. W., Marshall, J. A., Levenson, N. A., 

Elitzur, M., & Houck, J. R. 2007, ApJ, 655, L77 
Honig, S. F, Beckert, T, Ohnaka, K., & Weigelt, G. 2006, A&A, 452, 459 
Ivezic, Z., & Elitzur, M. 1997, MNRAS, 287, 799 
Ivezic, Z.et al. 2002, AJ, 124, 2364 

Ivezic, Z., Nenkova, M., & Elitzur, M. 1999, User manual for 
DUSTY, University of Kentucky Internal Report, accessible at 
http : // www . pa . uky . edu/~moshe/ dusty/ 

Jaffe, W. et al. 2004, Nature, 429, 47 

Kaspi, S., & Netzer, H. 1999, ApJ, 524, 71 

Kishimoto, M., Antonucci, R., & Blaes, O. 2005, MNRAS, 364, 640 
Krolik, J. H., & Begelman, M. C. 1988, ApJ, 329, 702 
Levenson, N. A., Sirocky, M. M., Hao, L., Spoon, H. W. W., Marshall, J. A., 
Elitzur, M., & Houck, J. R. 2007, ApJ, 654, L45 



Maiolino, R., Shemmer, O., Imanishi, M., Netzer, H, Oliva, E., Lutz, D., & 

Sturm, E. 2007, A&A, 468, 979 
Mason, R. E., Geballe, T. R., Packham, C, Levenson, N. A., Elitzur, M., 

Fisher, R. S., & Perlman, E. 2006, ApJ, 640, 612 
Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425 
Natta, A., & Panagia, N. 1984, ApJ, 287, 228 
Nenkova, M., Ivezic, Z., & Elitzur, M. 2002, ApJ, 570, L9 
Nenkova, M., Sirocky, M. M., Nikutta, R., Ivezic, Z., & Elitzur, M. 2008, 

ApJ, submitted 

Netzer, H. 1990, in Active Galactic Nuclei, ed. T. J. L. Courvoisier & 
M. Mayor, 57 

Ossenkopf, V., Henning, T, & Mathis, J. S. 1992, A&A, 261, 567 
Pier, E. A., & Krolik, J. H. 1992, ApJ, 401, 99 
Poncelet, A., Perrin, G., & Sol, H. 2006, A&A, 450, 483 
Rowan-Robinson, M. 1995, MNRAS, 272, 737 

Schartmann, M., Meisenheimer, K, Camenzind, M., Wolf, S., & Henning, T. 

2005, A&A, 437, 861 
Silva, L., Maiolino, R„ & Granato, G. L. 2004, MNRAS, 355, 973 
Sirocky, M. M., Levenson, N. A., Elitzur, M., Spoon, H. W. W., & Armus, L. 

2008, ApJ, 678, 729 
Spoon, H. W. W., Marshall, J. A., Houck, J. R., Elitzur, M., Hao, L., Armus, 

L., Brandl, B. R., & Charmandaris, V. 2007, ApJ, 654, L49 
Tomita, H. et al. 2006, ApJ, 652, L13 
Tristram, K. R. W. et al. 2007, A&A, 474, 837 
Urry, C. M., & Padovani, P. 1995, PASP, 107, 803 
Wolf, S. 2003, ApJ, 582, 859 



